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1. Introduction 

Lattice QCD is now established as a precision tool for calculations of the properties of 'gold- 
plated' hadrons. This enables us to test QCD at a much more stringent level than is possible for 
'QCD-inspired' models and is a necessary prerequisite to trust lattice QCD in more speculative 
calculations. It also enables us to determine accurately the parameters of the Standard Model 
related to quarks, including quark masses, the strong coupling constant and elements of the CKM 
matrix. 

Recently it has become possible to handle charm quarks accurately in lattice QCD using the 
Highly Improved Staggered Quark (HISQ) action [|l|] and this adds an additional new dimension 
to this programme. Since charm quarks straddle the region between light quark physics and heavy 
quark physics, special care must be taken. With the HISQ action charm quarks are treated in the 
same way as light quarks and this has a number of advantages, discussed below. Our recent results 
on charm-light decay constants [Q] show the power of this approach with results of comparable 
accuracy to those from light-light decay constants. Further results in charm physics are described 
in Section 2 and then a new method for determining the charm quark mass (to 1%) using our results. 
Preliminary results on the determination of ntb using the same method along with numerical results 
with the NRQCD action for b quarks are also given. An update on the determination of a s from 
Wilson loops completes section 3. The conclusions includes a summary plot of the current status 
of the gold-plated meson spectrum from lattice QCD. All of these results use Bayesian methods to 
allow a unified treatment of systematic errors from unknown higher order terms in fitting functions. 
Ignoring these terms clearly underestimates their effect; including them enables much more robust 
extrapolations and error estimates. How to do this is described in the Appendix. 

2. Charm physics with HISQ quarks 

Charm quarks are heavy, but not very heavy, and this makes a lattice QCD approach to them 
difficult. A heavy quark treatment effectively removes the mass as a dynamical scale, and thereby 
the discretisation errors linked to it. A light quark treatment allows the determination of the mass 
from the energy at zero momentum and, given an action with enough chiral symmetry, conserved 
currents that do not need renormalisation. As lattices become finer and finer, the charm quark mass 
in lattice units becomes smaller and the advantages of treating charm quarks as light relativistic 
quarks become more apparent. For example, on the MILC superfine ensemble (a fx 0.06fm) m c a 
is around 0.3, which is clearly much less than 1. The key issue is that of discretisation errors. 
These will appear as powers of m c a and can be extrapolated away, given results at enough values 
of a. However, if they are large i.e. low powers of m c a are present, so that a large extrapolation 
is necessary, they will cause significant errors in the extrapolated result. For a value of m c a = 0.4, 
tree-level errors at (m c a) 2 could cause a 20% error and a s (m c a) 2 a 6% error. It is therefore very 
important to use a highly improved action, in which these terms are removed, for an accurate 
continuum extrapolation. 

The HISQ action has been developed with exactly these issues in mind. It includes a further 
application of the 'FatV type gluon-link smearing beyond that used in the improved staggered 
(asqtad) action. As for that action it includes the improvement of the derivative using the Naik 3- 
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Figure 1: The difference in mass between the goldstone pseudoscalar and the next lightest taste, which is 
the local nongoldstone generated from the local temporal axial current for ss (rj s ) and cc (rj c ) and a variety 
of staggered actions. HISQ is the Highly Improved Staggered Quark action and HISQ2 and HISQ3 include 
one or two further levels of Fat7 smearing and associated reunitarisation and adjustment of the Lepage 
discretisation correction. 



link term which removes tree-level a 2 errors. For charm quarks the coefficient of the Naik term is 
further tuned nonperturbatively in HISQ to give 1 as the value for the 'speed of light'. This removes 
the leading (in terms of the velocity of the charm quark inside the meson) terms at a s {m c a) 2 and 
(m c a) 4 . Further details are given in 

The HISQ action is a staggered action and so has multiple 'tastes' of quarks and mesons made 
from them, whose differences vanish as a — > 0. One test of potential discretisation errors from 
'taste-changing' effects is to study the mass splittings between different tastes of pseudoscalar 
meson. These effects show up there and are very small elsewhere in the spectrum. Figure |] shows 
the splitting between the 'goldstone' t] c and T] s and the next lightest meson (which is made with 
the local temporal axial current) as a function of lattice spacing. The taste-splitting is smaller 
for charmonium than for strange-onium because the mass is heavier and both are clearly falling 
rapidly with a 2 . We also show the taste-splittings with actions called 'HISQ2' and 'HISQ3' in 
which further applications of the Fat7 smearing are made. It is clear that the taste-splittings can be 
reduced to negligible levels, even on quite coarse lattices, with this procedure. 

The HISQ action has enough chiral symmetry for conserved vector and partially conserved 
axial vector currents so that matrix elements can be calculated on the lattice with no need for renor- 
malisation. This removes a major source of systematic error. The HISQ action is also numerically 
very fast and so accurate results can be obtained quickly using variance reduction techniques such 
as the 'random wall' previously used for light quarks [[$J]. We fit correlators using a Bayesian multi- 
exponential approach |Q] taking account of oscillating states for the charm-light case (they are not 
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present for charmonium). The appendix describes the Bayesian approach in general terms. The 
form for the fitting function including normal (even) and oscillating (odd) states is 

G(t) = Y,c n (-l) nt (e- M "' +e- M "^) (2.1) 

We generally take our ground state energies and amplitudes from fits that include 4 normal (and 
4 oscillating if appropriate) exponentials. By this stage the fit results and errors are stable as a 
function of the number of exponentials, so it is irrelevant if more are added. 

We have tuned the charm quark mass in the HISQ action using the goldstone r] c mass on 
very coarse, coarse and fine MILC ensembles including 2+1 flavors of sea quarks with the asqtad 
action and various light quark masses. The scale is determined using MILC r\ja values [fij] and 
taking r\ = 0.321 fm [ph. We will also include here in some plots preliminary results from the 
MILC superfine (0.0036/0.018) ensemble. By combining HISQ charm quark propagators with 
light valence HISQ quarks we obtained accurate results for the mass of the D and D s mesons in [Q]. 
These masses have no free parameters because m c was already tuned from the r\ c and they provide 
a very strong test that charmonium and charm-light physics is handled simultaneously correctly in 
lattice QCD. Thinking of charm quarks as heavy quarks makes this a non-trivial test, not met by 
any continuum model of these systems. 

Our results for the D and D s decay constants, made in 2007, have caused a lot of interest 
this year following new results from CLEO-c ^ We calculate the decay constant from the 
matrix element of the partially conserved axial current in the standard way used for f n and fx- 
Discretisation errors are somewhat larger for D / D s than for Tt/K but clearly still under good control 
and extrapolation to the continuum limit from 3 values of the lattice spacing gives a 2% accurate 
result. Again we use Bayesian methods for a combined chiral and continuum extrapolation as 
described in the appendix. For D mesons the chiral extrapolation is significant; for D s mesons 
there is very little physical sea quark mass dependence, as expected. The decay constant is the 
quantity calculable in QCD that parameterises the annihilation rate of charged pseudoscalars to W 
bosons and thereby to leptons. Experimentalists can then also determine decay constants from the 
decay rate to leptons given a value for the appropriate CKM element from elsewhere. CLEO-c 
do this for the D and D s mesons (measuring the rate either to /xv or to TV and applying a small 
correction for electromagnetic effects) using V cs = Vud and V cc i = V us . Their new value for fo agrees 
well with our predicted result, their value for fu s disagrees at the level of 3 a where a is dominated 
by the experimental error coming from limited statistics. 

Figure || shows a compilation of lattice results for the decay constants, including an update 
from FNAL/MILC and a new value from n/ = 2 calculations by ETMC [@|. Also included 
are the new results from CLEO-c this year []6[ [7J, and older results on fj) s from BaBar [ fTO] ] and 



Belle [11]. There is clearly tension between the lattice results for fp s and the experimental ones. 
This is the first time that lattice QCD has disagreed with experiment on a gold-plated quantity, and 
over 15 such quantities have now been calculated accurately, for examples see figure [7|. It could be 



a harbinger of new physics [12]. At the very least it requires everyone to check their results and 
errors thoroughly. 

To this end, as well as running at a fourth superfine lattice spacing, we have been examining 
other gold-plated quantities in charm physics calculable using the HISQ action. Figure || shows a 
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Figure 2: A comparison of lattice results for the D and D s decay constant and experimental results obtained 
from the leptonic decay rate using CKM elements V CJ and V cc i from elsewhere. The FNAL/MILC results 
have been updated this year by [||] and the ETMC results are new and described in [^|. They include only 
2 flavours of sea quarks, so are not directly comparable to the results above. There is agreement between 
lattice and experiment for fo, but not for fn s . 



calculation of the decay constant of two vector mesons, the \ff and the . The decay constant can 
be extracted from the experimental leptonic width using 



— a 2 p 2 
3 u ^my 



tt 2 OE D e 2 o— (2-2) 



We compare results using different tastes of vector meson, 'local' and (taste-singlet) '1-link'. Nei- 
ther is the conserved current and so the lattice results had to be renormalised to compare to ex- 
periment. For the charmonium case this was done 'nonperturbatively' (i.e. using only continuum 
perturbation theory and not lattice perturbation theory) making use of current-current correlators as 
described in the next section. For the the renormalisation was done using 1-loop lattice perturba- 
tion theory, which apparently has small coefficients in the HISQ case. Agreement with experiment 
is clear, work is ongoing on the error budget. 



3. Determination of m c , rrif, and a s 

The accurate determination of quark masses is important for several continuum QCD calcula- 
tions. This is particularly true of the b and c quark masses whose uncertainty strongly affects the 
determination, for example, of V u }, from inclusive B — > % decays. The limitation on this determi- 
nation from standard lattice QCD methods is often the matching from the lattice bare mass to a 
continuum scheme such as MS. If lattice perturbation theory is used, a 2-loop determination of the 
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Figure 3: Results for the decay constant for the vector cc (J/y) and ss (0) using currents of different 
taste. The y/ case is renormalised 'nonperturbatively' using results from the comparison of the correlators to 
continuum perturbation theory. The <j> case is renormalised using one-loop lattice QCD perturbation theory 
(no error is included for unknown higher orders in the renormalisation in this plot because such an error is 
correlated between the points at different lattice spacing). No dependence on the light sea quark mass is seen 
in either of these quantities, as expected. 



matching factor must be done and this is hard. Nevertheless encouraging results can be obtained in 



the case of the charm quark using the HISQ action Q13|]. 

Here we describe a new method which takes a very different approach. It requires a com- 
parison of the continuum extrapolation of zero-momentum lattice charmonium correlators to high 
order continuum QCD perturbation theory, and the work was done in collaboration with Chetyrkin, 



Kiihn, Steinhauser and Sturm who performed the continuum calculations [14]. The comparison is 
done through the nth 'time-moments' of the correlators (defined below) which can be related to nth 
derivatives with respect to qo, evaluated at q 2 = 0, of the polarisation function of an external current 
coupled to a heavy quark loop. This latter quantity is calculable in continuum perturbation theory, 
provided that n is not too large. Our results are most accurate for the goldstone pseudoscalar (r/ c ) 
correlator. We can simply multiply by the square of the bare charm mass to define an ultra-violet 
finite unrenormalised (because of the PCAC relation) current-current correlator: 

G(t) = a 6 £(am 0:C ) 2 < 0\j 5 {x,t)j s {0,0)\0 > (3.1) 

X 

and calculate time-moments as: 

G n = J £(t/a) n G(t). (3.2) 
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Here t goes from -T /2 to +T /2 on the lattice. The comparison to the continuum then becomes 

G ' i(a - 0) - (am^))-* • (3-3) 

where g„ is known through third-order in a s for low values of n (4, 6 and 8). The approach to the 
continuum limit is improved if tree level discretisation errors are removed. This is done by dividing 



by the tree level result on both sides of eq |3.3| . The tree level lattice result is simply obtained by 
calculating in the free case. Tuning errors in m c and scale setting errors are also ameliorated by 
multiplying by factors of the lattice r] c mass. Then the continuum extrapolation is actually done for 

_ am nc G n l/(n _ 4) 

K "-~ (0) ^ (0)' 
2am pole,c °n 

for n > 6. For n = 4 there are no factors of m c in eq. [O], so we cannot use this to obtain m c , but 
a s can be determined as described below. Extrapolation of the reduced time-moments, R n , to the 
continuum limit again uses Bayesian methods (see appendix): 

R n (a) =R n (0)(l +cia s (am c ) 2 + C2<x s (am c ) 4 + CT,a s (am c ) 6 -\ (3.5) 

We also include dependence on the sea quark masses which has negligible effect and is omitted 
here for clarity. Having 4 values of the lattice spacing from 0.06fm (MILC superfine) to 0.1 5fm 
(MILC very coarse) enables a very accurate continuum result to be obtained. The ratio of m c to the 
experimental r\ c mass is determined for different values of n using: 

r n (a m ,n/m c ) 

K " = o — 777w (.3-0) 



where r n is g„ from eq. ^3J divided by the continuum tree level result. A 1 % final error is obtained 



by averaging results from n = 6 and n = 8. The error is dominated by uncertainties in the scale 



setting and in the perturbation theory. Further details are given in [|14|], which includes also determi- 
nations of m c from temporal axial and vector currents and a full error budget. Figure |, from [0], 
summarises the results. Our final number is m c (3GeV) = 0.986(10)GeV. This method has been 



previously applied [15] to determine m c to 1% using experimental data on R(e + e~ — > hadrons). 
The continuum result is m c (3GeV) = 0.986(13)GeV. The fact that the lattice and continuum de- 
terminations of m c agree at this level of precision is a strong statement about how well we can 
handle charm quarks in lattice QCD using the HISQ action. 

Another application of these techniques is for the 'non-perturbative' determination of renor- 
malisation factors for the cases where a nonconserved current is used in the correlator. Different 
moments have the same Z factors, but different powers of m c appear. It is therefore possible, by 
taking ratios of adjacent moments to appropriate powers, to isolate Z. These factors were used in 
the calculation of the leptonic width of the y above. 

Since the 4th moment has no powers of m c it is possible to relate its continuum extrapolation 
directly to a perturbative expression and use this to determine a s . We can do this also, but less 
accurately, for the ratios of powers of the 6th and 8th moments chosen to cancel powers of m c . 
The value we obtain in the MS scheme is a s (Mz) = 0.1174(12), in good agreement with other 
determinations. More details are given in JT4|]. 
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Figure 4: m c (ji) for /I = 3 GeV and rif = 4 flavours, from different moments of correlators built from 
four different lattice operators. Top left is the local (goldstone) pseudoscalar, top right the local temporal 
axial vector and bottom left and right are the ' 1-link' vector and local vector respectively. The grey band is 
0.986(10) GeV, which comes from the first two moments of the local pseudoscalar. 
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Figure 5: Results for m& is the MS scheme at its own scale obtained from moments of vector (T) current- 
current correlators in lattice NRQCD matched to continuum perturbation theory as described in the text. The 
shaded band is 4.18(4) GeV. At small n the errors are dominated by relativistic and discretisation errors and 
at large n by nonperturbative uncertainties. 

The method of moments of current-current correlators can also be applied to extract mi, from 
T and r\b correlators using NRQCD b quarks on the MILC configurations []T(J. Now there is no 
conserved current in the annihilation channel and so ratios of moments must be taken to eliminate 
the factors of Z. Moments for small values of n suffer from discretisation and relativistic correc- 
tions, but the region of validity of the perturbative calculation extends to larger values of n. A tree 
level analysis of lattice results, compared to tree level full QCD, is useful in understanding the 
systematic errors. Figure || shows preliminary results for from T correlators. It seems likely 
that we can achieve a 1-2% result for mt from this method. Work is ongoing. 

Because the HISQ action is an appropriate action for light quarks too, we can take advantage 
of cancellations of systematic errors in ratios of charm quantities to light ones. This allows us to 
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lattice strange quark mass using the HISQ action and variants on it. HISQ2 (double hisq) and HISQ3 (treble 
hisq) are described in the text. Matching factors to the continuum MS masses cancel in this ratio and allow 
an accurate m s to be determined. 



leverage accuracy in light quark quantities from our charm results (and is therefore, it turns out, 
a very useful alternative to the more traditional approach of using the same action for charm and 
bottom and leveraging accuracy in bottom physics from charm [^]) 

Figure ^ shows the ratio of the charm quark mass to the strange, obtained with HISQ, HISQ2 
and HISQ3 actions on the MILC configurations. The matching to the MS scheme cancels in this 
ratio, when extrapolated to a = 0. We are in the process of using this for an alternative accurate 
determination of m s . There are clear discretisation errors in the ratio in the HISQ case and these 
are significantly reduced with the more highly smeared actions. Again work on this is ongoing. 

Finally, we give an update to the determination of a s from lattice QCD using the perturba- 
tive expansion of 22 small Wilson loops or loop ratios JTv|] . This is improved from our previous 
determination [18] in a number of ways and obtains a result which is slightly more accurate and 
la higher at 0.1184(9). Our new result includes determination from 5 values of the lattice spac- 
ing (adding the MILC superfine and very coarse) and makes use of the MILC r\/a values to fix 
the ratios of lattice spacings between different ensembles more accurately. We have also fitted the 
nonperturbative chiral corrections, ie. the light sea quark mass dependence of the logs of the Wil- 
son loops. This is a small effect but certainly required by the data when using results at multiple 
sea quark masses and lattice spacing values. The formula for the log of a small Wilson loop then 
becomes 

log(W) =Y J c n av(d/a)(l + c m a(2m l +m s ) + ■■■). (3.7) 

n,m 

The measured ensemble average value on the MILC ensembles is inserted on the left and the equa- 
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tion inverted to obtain a,y(d/a). The scale for a, d/a is set using the BLM prescription, modified 
where necessary [[l9|]. c\, C2 and C3 are known from the numerical evaluation of Feynman dia- 
grams in lattice QCD perturbation theory, c\ ...c\o are constrained from our Bayesian analysis (see 
appendix) which fits results for a particular log(W) from the 5 different lattice spacings using a 
common ay, which runs perturbatively between the different scales. We are able to constrain the 
higher order terms somewhat more accurately than in our previous calculation because of having 
the additional superfine results, and these are largely responsible for the shift upwards of our num- 
ber. The inclusion of unknown higher order terms is necessary for an accurate analysis, and the 
Bayesian approach allows us both to parameterise our uncertainty about them as well as allowing 
the data to constrain them where it can. They give sizeable corrections and without them a poor fit 
would be obtained across the multiple lattice spacings. The Bayesian fit also allows us to include 
uncertainties in the scale (both statistical and systematic) and analyse the effect of gluon conden- 
sates of various kinds. Further details are given in [|l7|]. See also [20] for a somewhat different 
analysis of the perturbation theory for a subset of the Wilson loops and which gives a result in good 
agreement. 

Note the difference in methodology between the two determinations of a s that we give here. 
The first (moments of current-current correlators) uses a quantity which is defined in the contin- 
uum and has lattice artefact (discretisation) errors which must be extrapolated away before being 
compared with continuum perturbation theory. The second (small Wilson loops) uses a quantity 
defined on the lattice and calculated in lattice perturbation theory which includes all lattice arte- 
facts, order by order in ay. Discretisation errors only enter here through the presence of such errors 
in the quantity used for scale setting in the scale for ay . The second method then has an advantage 
provided that lattice perturbation theory can be done to high enough order, and our errors reflect 
this. 



4. Conclusions 



Precision lattice QCD calculations continue to produce valuable results, building on pl[], and 
we give new examples here. With the advent of the HISQ action, charm physics has become an 
excellent testing ground for lattice QCD and QCD. The accuracy of b physics now needs to be 
improved to the same level and work on this is ongoing [22]. Figure [7] shows the current status of 
the gold-plated meson spectrum from HPQCD calculations. 
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Appendix : Constrained Fits and Error Budgets 

Constrained fitting, with Bayesian priors, is the most reliable tool for analyzing systematic 
errors associated with correlator fits and continuum, chiral and other extrapolations that involve 
a large or infinite number of increasingly unimportant terms [p3p. Fitting lattice results to a chi- 
ral expansion truncated at next-to-leading order, to take a widely used example, provides little 
unambiguous information about the potential impact of higher-order terms in the expansion. In a 
constrained fit, an arbitrary number of higher-order terms can be included and their potential impact 
on systematic errors easily quantified. 

The key to constrained fitting lies in the Bayesian priors that are included in the % 2 function 
that is minimized in the fit. Typically one is trying to fit a set of data points, say F,- ± Oy i at points 
X{ for i = 1 . . . Ny (for example, simulation results for correlators as a function of t or values for 



R n defined in eq. [3^ or log(VK) for different lattice spacings), to a function y{x\c) that depends 
upon a set of fit parameters cj for j = I ...N c (for example, amplitudes and masses in the sum 
of exponentials in eq. 2.1 or coefficients of the powers of (am c ) 2 in eq. or of ay{d/a) or 



a(2mi +m s ) in eq. 3.7). A complication arises when in principle N c is infinite, but in practice 
only the first few cys contribute appreciably. We usually can estimate roughly how many terms are 
needed, because we have prior knowledge about the order of magnitude of the CjS. The challenge 
is to incorporate systematically this prior knowledge into the fitting process so that the impact of 
unimportant or marginally important terms can be reliably assessed and quantified. 

This information enters a constrained fit through the Bayesian priors for the fit parameters. In 
such a fit, we vary the fit parameters cj to minimize an augmented x 2 function, 

Z 2 (c)=L (y< " y( f ;C))2 +Ig^, (4-D 

i °Y, j 

where there is one prior SXc for every fit parameter. The priors incorporate our prior knowledge 
about the fitting parameters. Typically one uses a Gaussian prior, 

SXI = ^P^, (4-2) 



to constrain the fit parameter cj to the vicinity of cj ± a c . . 
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The mean cj and width a Cj encode our prior knowledge about fit parameter cj. The data being 
fit will have a lot to say about some parameters, and almost nothing to say about others. In the 
former case, where a parameter is very sensitive to the data, usually the fit result for that parameter 
has an uncertainty that is much smaller than the width of its prior, and both the mean and standard 
deviation that come out of the fit are almost independent of the mean and standard deviation that 
went into the prior. In the other extreme, where the data is largely insensitive to a parameter, the 
mean and standard deviation that come out of the fit are approximately equal to the mean and 
standard deviation put into the prior. In the first case, the data gives us new information beyond 
what we put into the priors (that is, beyond what we knew before doing the fit); in the second case, 
the data adds nothing to our knowledge about the parameter. Parameters of the second kind, which 
have little relevance to the data, are good examples of nuisance parameters in statistics. 

Bayesian analysis provides a logically coherent framework for this kind of analysis [23]. Its 
power lies in the fact that we can include an arbitrary number of parameters without destabilizing 
the fit. This is as it should be. If a set of parameters is irrelevant because their corresponding terms 
in y{x;c) are negligible, we should be free to include them or not in an analysis. It should make 
no difference. If they are marginally important it should have a marginal impact on the fit. This is 
precisely the situation when proper priors are included. 

As one adds cjs in the fitting function, the quality of the fit should initially improve (that is 
% 2 decrease), and usually the errors on fit results increase. Eventually, however, neither the fit nor 
the errors change when further parameters are added. This is the point at which the parameters 
become insensitive to the data (usually because the data is insufficiently accurate to resolve them). 
This typically occurs around the point where % 2 /Ny falls below one. It is important to add terms 
up to the point where further terms don't change the results, so that systematic errors are not 
underestimated. Terms can be added beyond this point but there is little merit in this since they 
have no effect on fit results (means or errors). 

Means and standard deviations for the fit parameters, and functions of them, are obtained from 

i 

'g 

linear in the variances that appear in the various terms in the % 2 function: 

ot™L d ri a T, + L d ctfr (4 " 3) 

i j 

The first sum, for example, is the contribution to the error in g(c) coming from statistical errors 
in the data. The second sum is the contribution from uncertainties in the fit function. This sort 
of information is obviously very useful when thinking about improvements to an analysis — for 
example, in deciding how much improvement would come from better statistics for the data. 

To isolate the part of the total error a g that comes from, for example, the statistical uncertainty 
in all the YjS, the fit is rerun but with the corresponding variances in the % 2 function rescaled by a 
factor / close to one (/= 1.1 or 1.01, for example): 

ol^fo 2 (4.4) 

fori = I... N r . Then 

o g (f) 2 -o 2 (f-- 
/-l 



the minimum x 2 in the usual fashion. It is often useful to decompose the error o g for for some fit 
result g(c) into component parts. When errors are small (as here), the variance a 2 is approximately 



£^<4 (4-5) 
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The square root of this quantity is the part of the total error that comes from statistical uncertainties 
in the Yj. This procedure can be repeated for each prior or group of priors in the x 2 function, 
thereby generating a complete error budget for the g(c). The sum of the variances obtained in this 
way for each part of the total error should equal cr?; if it does not, errors may not be small enough 



to justify the linear approximation in 4.3 



Most of the uncertainties in a standard lattice analysis can be pushed into the constrained fit, to 
facilitate consistent treatment of all systematic and statistical errors. For example, the independent 
variables X; might have errors: X,- = X, ± Ox r The determination of the lattice spacing will typically 
have a systematic error that grows with the lattice spacing that can be taken into account this way, 
for example. We do this in the analysis of Wilson loops to derive ay since no continuum limit 
is taken there but the systematic errors in the comparison of scales between ensembles is very 
important. Another example is that of chiral extrapolation using quark masses, where the quark 
masses have to be run to a common scale and will therefore have systematic errors associated with 
them. In these cases we treat each X, as another fit parameter, to be varied in the fit, and include the 
following priors in the % 2 function: 

8X 2 X = £ (4-6) 
i a x, 

One of the outputs from the fit might be improved values for the X,. 

There are usually a variety of other parameters that introduce errors into a lattice analysis. 
For example, we take the lattice spacing from MILC data for r\/a and use our T determination of 
r\ [Q, |16| ]. Each of these has errors. The {r\/a)i for each data set and the overall r\ are all treated as 
fit variables when we fit our lattice data to a function that depends on the lattice spacing, and each 
has a prior: 



Here {r\jd) i and C7( ri / a ). are taken from MILC results, while r\ ± a n = 0.321(5) fm. 

Any quantity that enters the analysis and has uncertainties is treated as a fit parameter, with a 
prior that incorporates whatever is known about that quantity. In this way all uncertainties can be 
analyzed simultaneously within a single constrained fit. 



Occasionally the difference in i.5 comes out negative. This could be caused by instabilities in the fit, in which case 
changing / might change the sign. It is possible, however, for one of the coefficients in 13 to be negative. In such cases 
we use the absolute value of 4.5 for the variance. 
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